* Settings
version 16
do "$SSDIMed/scripts/_auxiliary/_project_settings.do"


* -------------------------------------------------------------------------------------------------
* Prepare data for figures
* -------------------------------------------------------------------------------------------------

* DI entrant outcomes, by x, where x = 
*   i.covstart_year          = calendar year of Medicare entry
*   c.unemprateatapp         = unemployment rate at time of application
*   i.unemp_q20              = ventiles of unemprateatapp
*   i.age_year_covstart_fill = age in years at Medicare entry

* Sample
*   main: individuals who enter Medicare at ages 20-62, and annual Medicare spending after first (usually partial) year of Medicare coverage

* Sample
*   main: individuals who enter Medicare at ages 20-62, and annual Medicare spending after first (usually partial) year of Medicare coverage
*   2262: main, but individuals who enter Medicare at ages 22-62
local sample 2262

* Sample definitions:
* - main (all ages, all years)
local sample_main 1
* - ages 22-62
local sample_2262 inrange(age_year_covstart_fill, 22, 62)
* - ages 51-52
local sample_5152 inrange(age_year_covstart_fill, 51, 52)
* - first full year
local sample_ffy (rfrnc_yr == init_yr + 1)

* Which controls were included in each regression
* y = incidence
*   local controls_01 i.age_year_covstart_fill
*   local controls_02 i.fipscounty_init_g##i.age_year_covstart_fill

* DI incidence
local y entry
foreach x in age_year_covstart_fill {  
  * Load data
  use "$SSDIMed/data/analysis/county-month-age_entry_sample-main.dta", clear
  gisid fipscounty_firstnm covstart_month age_year_covstart_fill, missok
  
  * Confirm original sample size
  sum count
  assert r(sum) == 16240776
  qui sum age_year_covstart_fill
  assert r(min)==20 & r(max)==62
  
  * Select sample
  assert `sample_`sample'' == inrange(age_year_covstart_fill, 22, 62)
  keep if `sample_`sample''
  sum count
  assert r(sum) == 15790262

  * Collapse
  gcollapse (mean) nat_unemp_mean=unemp_rate_national_atapp  incidence=incidence_pop_age_atapp [aw = pop_age_atapp], by(`x') fast
  
  * Save
  tempfile `x'_`y'
  save ``x'_`y'', replace
  list, sep(0)
}


* DI incidence vs. unemployment
local y incidence_pop_age_atapp
local x c.unemp_rate_county_atapp
local clust county_mofd
if 1 {
  * Control variables
  local ctrl 02
  
  * Confirm sample
  assert "`sample'" == "2262"
  
  * Load estimates of cyclical entry by age
  estimates use "$SSDIMed/results/estimates/x-`x'/x-`x'-age_year_covstart_fill_y-`y'_controls-`ctrl'_cluster-`clust'_sample-`sample'.ster"
  regsave, ci detail(all)
  
  * which observations and variables to keep
  drop if var == "_cons"
  keep var coef stderr ci_lower ci_upper
  
  * numeric level for the age
  gen age = substr(var, 1, strpos(var, "."))
  replace age = regexr(age, "(b\.)|(\.)", "")
  destring age, replace
  order age, after(var)
  
  * Merge in incidence levels
  gen age_year_covstart_fill = age
  merge 1:1 age_year_covstart_fill using `age_year_covstart_fill_entry', assert(match) nogen
  
  * Incidence response as % of incidence level
  foreach var in coef stderr ci_lower ci_upper {
    replace `var' = (`var'/incidence) * 100
  }
  
  list, sep(0)
}


* ---------------------------------------------------------------------------------------
* Graph settings
* ---------------------------------------------------------------------------------------

* Set graph style to project default settings
*	Pass as arguments: height width (default 3.5in 6.5in)
do "$SSDIMed/scripts/_auxiliary/_project_grstyle.do" 3.1in 6.5in

local red "179 0 0"
local green "75 115 47"
local blue "48 84 150"
local bluegray "126 153 180"
local brown "140 045 004"
local orange "237 125 49"
local tan "210 180 140"

* Entry
local p 2
grstyle set color "`brown'": p`p'
grstyle set symbol i: p`p'
grstyle set symbolsize 2.2pt: p`p'
grstyle set lpattern solid: p`p'
grstyle set linewidth 1.8pt: p`p'

* Arrows and lines
local p 5
grstyle set color gs8: p`p'
grstyle set lpattern solid: p`p'
grstyle set linewidth vvthin: p`p'


* -------------------------------------------------------------------------------------------------
* Figure: Cyclicality of DI entry by age at Medicare entry
* -------------------------------------------------------------------------------------------------

* Describbe data
list, sep(0)
assert inrange(age, 22, 62)

* xvar is age at entitlement, taken to be two years before Medicare coverage starts
cap drop xvar
gen double xvar = age - 2
assert inrange(floor(xvar), 22 - 2, 62 - 2)
qui glevelsof xvar
assert r(J) == (62 - 22 + 1)


if 1 {

* Make min and max labels align on both y axes
* Stata tip 93: Handling multiple y axes on twoway graphs: https://journals.sagepub.com/doi/pdf/10.1177/1536867X1101000411

* How many major ticks?
* How much should range extend beyond ticks (as fraction of tick step size)
local y_ticks 8
local y_margin_lo 0
local y_margin_hi 0.4

* y1 tick lo, step size, format
local y1_tick_lo 0
local y1_step 1
local y1_fmt "%9.2g"
local y1_unit "%"

* Construct y labels and scales from parameters above
forvalues y = 1/1 {
  local y`y'_tick_hi  = `y`y'_tick_lo' + `y`y'_step' * (`y_ticks' - 1)
  local y`y'_range_lo = `y`y'_tick_lo' - `y`y'_step' * `y_margin_lo'
  local y`y'_range_hi = `y`y'_tick_hi' + `y`y'_step' * `y_margin_hi'
  di "y`y'_tick_hi:  `y`y'_tick_hi'"
  di "y`y'_range_lo: `y`y'_range_lo'"
  di "y`y'_range_hi: `y`y'_range_hi'"
  
  local y`y'_label ylabel(`y`y'_tick_lo'(`y`y'_step')`y`y'_tick_hi', axis(`y') format("`y`y'_fmt'") noticks nolabels labgap(1pt)) 
  local y`y'_scale yscale(`yalt' range(`y`y'_range_lo' `y`y'_range_hi') axis(`y') noline) 
  di `"y`y'_label:  `y`y'_label'"'
  di `"y`y'_scale:  `y`y'_scale'"'
}

* x-axis settings
local xlo 18
local xhi 61
local xscale xscale(range(`xlo' `xhi'))
local x_axis xlabel(20(5)60) `xscale' xtitle("Age at entry", margin(t=+.5))

* style and color for y-axis title/labels
local y1_title_style it
local y1_title_color gs8
local y1_label_style it
local y1_label_color gs8
local y1_label_size 9pt

local text_style sf
local text_color black

* y-axis labels placed above grid rule (build manually)
local y1_x `xlo'
local y1_place ne
local y1_just left

local y 1
local y`y'_la
foreach tick of numlist  `y`y'_tick_lo'(`y`y'_step')`y`y'_tick_hi' {
  if `tick' == `y`y'_tick_hi' {
    * white textbox to go behind label
    local y`y'_la `y`y'_la' text(`=`tick' + `y`y'_step'/100' `y`y'_x' " ", yaxis(`y') place(`y`y'_place') just(`y`y'_just') size(`y`y'_label_size') box lwidth(none) color(`y`y'_label_color') width(5) height(3) bcolor(white%70))
    *label
    local y`y'_la `y`y'_la' text(`tick' `y`y'_x' "{`y`y'_label_style':`=string(`tick', "`y`y'_fmt'")'`y`y'_unit'}", yaxis(`y') place(`y`y'_place') just(`y`y'_just') size(`y`y'_label_size') box lwidth(none) color(`y`y'_label_color') bcolor(none))
  }
  else {
    local y`y'_la `y`y'_la' text(`tick' `y`y'_x' "{`y`y'_label_style':`=string(`tick', "`y`y'_fmt'")'}", yaxis(`y') place(`y`y'_place') just(`y`y'_just') size(`y`y'_label_size') box lwidth(none) color(`y`y'_label_color') bcolor(none))
  }
}

* y-axis titles, horizontal orientation (build manually)
local y 1
local y1title text(`=`y`y'_tick_lo' + `y`y'_step'*(`y_ticks' - 1 + `y_margin_hi')' `xlo' "{`y`y'_title_style':Entry effect as a percentage of total age-specific entry}"
local y1title `y`y'title', yaxis(`y') place(`y`y'_place') just(`y`y'_just') box lwidth(none) color(`y`y'_title_color') bcolor(white) margin(b=+.5))
local y1title `y`y'title' ytitle(" ", axis(`y') margin(l=-8))

* y-axis settings
local y_axis1 `y1_label' `y1_la' `y1_scale' `y1title'

* Plot lines
local line_count (connected coef xvar, yaxis(1) pstyle(p2))
local ci_area (rarea ci_lower ci_upper xvar, sort pstyle(p2area))

* Label line 1
* local age 51
* sum coef if age == `age'
* local coef = r(mean)
* sum stderr if age == `age'
* local stderr = r(mean)
* local text_1 text(1.22 `=`age'-1.6' "Age `age' percentage effect: `=string(`coef', "%9.1f")' (se = `=string(`stderr', "%9.1f")')"
* local pci_1 (pci 1.2 `=`age'-1.5' `coef' `=`age'', pstyle(p5))

* Label line 2
* local age 52
* sum coef if age == `age'
* local coef = r(mean)
* sum stderr if age == `age'
* local stderr = r(mean)
*local text_2 text(6.22 `=`age'-1.6' "Age `age' percentage effect: `=string(`coef', "%9.1f")' (se = `=string(`stderr', "%9.1f")')"
*local pci_2 (pci 6.2 `=`age'-1.5' `coef' `=`age'', pstyle(p5))

* Eligibility changes
local text_grid text(`=`y1_tick_lo' + 0.58*`y1_step'' 52.5 "Eligibility" "relaxes", place(c) just(center) box lwidth(none) bcolor(white%94))
local xlines_grid xline(50) xline(55)
local pci_grid (pci 0.6 50.6 0.59 50, pstyle(p5)) (pci 0.6 54.4 0.59 55, pstyle(p5))

local legend legend(off)
local region graphregion(margin(t=5 b=-1))

* twoway `ci_area' `line_count' `pci_grid' `pci_1' `pci_2', `x_axis' `y_axis1' `legend' `region' `text_grid' `xlines_grid' yline(0, lpattern(solid) lwidth(thin) lcolor(gs8)) `text_1' `text_2'
twoway `ci_area' `line_count' `pci_grid', `x_axis' `y_axis1' `legend' `region' `text_grid' `xlines_grid' yline(0, lpattern(solid) lwidth(thin) lcolor(gs8))

* Export graph
cap mkdir "$SSDIMed/results/figures"
graph export "$SSDIMed/results/figures/entry_unempatapp_ageatentry_pct-`sample'.pdf", as(pdf) fontface("Source Sans Pro") fontdir("$SSDIMed/data/raw/fonts") replace

}






** EOF
